"""
该程序为粒子群算法求解第二问（连续型非线性规划）
主要改进：
1 GBest PSO（Global Best Particle Swarm Optimization）
GBest PSO 是 PSO 的一种变体，其中包括全局最优粒子（Global Best Particle）
2 融入先验知识： 30% 的初始粒子人为指定 70% 粒子随机生成（可以帮助收敛）
3 根据历史信息限定决策变量范围（上下界）
"""

import numpy as np
import matplotlib.pyplot as plt
import pyswarms as ps
list_l = ['花叶类' ,'花菜类' ,'水生根茎类' ,'茄类' ,'辣椒类' ,'食用菌' ]

predit_buy =[
    [3.285864, 3.2921748, 3.2889733, 3.285188, 3.2851105, 3.2964268, 3.2876368] ,
    [7.7414317, 7.763459, 7.814592, 7.794937, 7.747068, 7.810813, 7.7633805] ,
    [12.018651, 11.912668, 12.027704, 11.941088, 11.92054, 12.118359, 11.972251] ,
    [4.5562034, 4.601929, 4.5483465, 4.549116, 4.532483, 4.539543, 4.601603] ,
    [3.7067149, 3.65774, 3.6644902, 3.6755412, 3.6658049, 3.6834998, 3.6471841] ,
    [4.015016, 4.075036, 4.025253, 4.0783653, 4.0397897, 4.0211616, 4.0569587] ]# 预测进价
predit_sale =[
    [190.75572, 189.94437 ,189.2342 , 190.34938 ,189.08669, 188.56415 ,190.10588],
    [28.618061, 28.872581, 28.873682 ,28.74203 , 28.776909, 28.86964 ,28.661997],
    [26.962397, 27.805391, 27.65219, 27.210875, 27.88252, 27.764929, 27.24346 ],
    [32.29913, 31.795496 ,32.425, 31.649815, 31.683603 ,31.381622, 31.967655],
    [102.14867, 101.67641, 102.01936, 102.55891 ,102.07538 ,102.73362 ,102.2766 ],
    [56.307552, 57.39569 , 58.130955, 56.80816, 56.4629, 57.186737, 57.953686] ]# 预测销量

predit_omega =[0.7]*6 # 折扣
predit_gama = [0.1283, 0.1551, 0.1365, 0.0668, 0.0924, 0.0945] # 损耗率

day = 6

ini_pos = [0.6, 0.4, 0.3, 0.6, 0.9, 0.6,150, 30, 20, 20, 100, 60] #决策变量：前六个是利润率，后六个是进货量
ini_pos = np.array(ini_pos)
n_particles = 1000 #粒子个数
n_dimensions = 12 #样本维度

#设置决策变量上下限便于收敛
lower_bound = np.array( [0.3, 0.2, 0.2, 0.3, 0.3, 0.3, 10, 0, 0, 10, 10, 20] )
upper_bound = np.array( [1.2, 0.98, 0.8, 0.9, 1.2, 1.0, 450, 90, 75, 60, 300, 250] )
bounds = (lower_bound, upper_bound)

weight_ini = 0.3 # 给定的初始值占总粒子的比例

#随机生成初始粒子
pos_given = np.random.uniform(
    low=lower_bound, high=upper_bound, size=(int(n_particles * weight_ini) , n_dimensions))
pos_given = 0.8 * pos_given + 0.2 * ini_pos
pos_given = np.clip(pos_given,lower_bound,upper_bound)
pos_random = np.random.uniform(
    low=lower_bound, high=upper_bound, size=(int(n_particles * (1-weight_ini)), n_dimensions))
Initial_pos = np.vstack((pos_given, pos_random))

#与粒子群无关的预测修正函数
def modify(y,x,idx):
    if idx == 0:
        sale = -1.91 * x + 1.03 * y + 3.82
        if sale < 0:
            return 0
        else:
            return sale
    if idx == 1:
        sale = -0.25 * x + 0.94 * y + 1.07
        if sale < 0:
            return 0
        else :
            return sale
    if idx == 2:
        sale = -0.03 * x + 0.97 * y + -0.98
        if sale < 0:
            return 0
        else :
            return sale
    if idx == 3:
        sale = -0.07 * x + 0.99 * y + 1.70
        if sale < 0:
            return 0
        else :
            return sale
    if idx == 4:
        sale = -0.029 * x + 0.99 * y + 2.54
        if sale < 0:
            return 0
        else :
            return sale
    if idx == 5:
        sale = -0.51 * x + 1.01 * y + 2.32
        if sale < 0:
            return 0
        else :
            return sale


# 目标函数，x是粒子
def Objective_function(x):
    """
    目标函数
    :param buy: 进价
    :param sale_lstm: 销售量（lstm 预测结果）
    :param alpha: 加成率（利润率）(变量)
    :param beta: 进货量 (变量)
    :param gama: 损耗率
    :param omega: 折扣
    :return: 收益
    """
    profit_list = []
    for i in range(n_particles) :
        x_new = x[i]
        profit = 0
        for idx in range(6):
            alpha = x_new[idx] # 利润率
            beta = x_new[idx + 6] # 进货量

            buy = predit_buy[idx][day] # 预测进价
            sale_lstm = predit_sale[idx][day] # 预测销量
            omega = predit_omega[idx] # 折扣
            gama = predit_gama[idx] # 损耗率

            sale_price_normal = buy * (1 + alpha) # 好货的售价
            sale_price_discount = buy * (1 + alpha) * omega # 差货的售价

            good = beta * (1 - gama) # 好的进货量
            bad = beta * gama # 差的进货量

            sale_modify = modify(sale_lstm, sale_price_normal,idx) # 预测的销量

            w1 = sale_modify * (1 - gama) * sale_price_normal + sale_modify * gama * sale_price_discount
            w2 = beta * buy

            #惩罚函数，惩罚进货量太大的情况
            if beta <= sale_modify:
                profit += (w1-w2) - (sale_modify - beta) * 20
            else:
                profit += (w1 - w2)

        profit_list.append(-profit)

    #优化目标是让profit最大
    return profit_list


options = {'c1': 0.5, 'c2': 0.5, 'w': 0.6} # 个人 社会 继承
optimizer = ps.single.GlobalBestPSO(n_particles=n_particles, dimensions=n_dimensions, options=options , bounds=bounds,init_pos=Initial_pos)

#best_cost是决策变量的取值，best_position是优化目标的取值，iters是迭代次数，verbose为true代表输出决策变量取值
best_position, best_cost = optimizer.optimize(Objective_function, iters=300, verbose= True)

#后续指标计算，与粒子群无关
def pre_sale_mount (x1,x2):
    # x1 是销售价格
    # x2 是 LSTM 预测的销量
    list_sale = []
    for idx in range(6):
        x = x1[idx]
        y = x2[idx]
        sale = modify(y,x,idx)
        list_sale.append(sale)

    return list_sale

def count_how_much(x,y):
# x 是利润率
# y 是进价
    list_how_much = []
    for i in range(6):
        list_how_much.append((x[i] +1) * y[i])

    return list_how_much
print('-------------------------------------')
print(" 最大利润：", -(best_position))

print(" 进货量（变量）：", best_cost[-6:].tolist())

print(" 预测进价：",[row[day] for row in predit_buy])
sale_price = count_how_much(best_cost[:6].tolist(),[row[day] for row in predit_buy] )
print(" 销售价格：",sale_price)
print(" 利润率（变量）：", best_cost[:6].tolist())

# print("\nLSTM 预测销量：",[row[day] for row in predit_sale] )
print(" 修正之后的销量：",pre_sale_mount(sale_price, [row[day] for row in predit_sale] ))